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Abstract 

When gravitational fields are at their strongest, the evolution of spacetime 
is thought to be highly erratic. Over the past decade debate has raged over 
whether this evolution can be classified as chaotic. The debate has centered 
on the homogeneous but anisotropic mixmaster universe. A definite resolu- 
tion has been lacking as the techniques used to study the mixmaster dynamics 
yield observer dependent answers. Here we resolve the conflict by using ob- 
server independent, fractal methods. We prove the mixmaster universe is 
chaotic by exposing the fractal strange repellor that characterizes the dynam- 
ics. The repellor is laid bare in both the 6-dimensional minisuperspace of the 
full Einstein equations, and in a 2-dimensional discretisation of the dynam- 
ics. The chaos is encoded in a special set of numbers that form the irrational 
Farey tree. We quantify the chaos by calculating the strange repellor's Lya- 
punov dimension, topological entropy and multifractal dimensions. As all of 
these quantities are coordinate, or gauge independent, there is no longer any 
ambiguity-the mixmaster universe is indeed chaotic. 
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Einstein's theory describes gravity from falling apples to the collapse of an entire universe. 
When space is most strongly deformed, Einstein's nonlinear theory may be fundamentally 
chaotic. Chaotic systems are complex but knowable. A paradigm shift from simplicity to 
complexity is needed. The attempt to track each individual universe is abandoned. Instead, 
as with thermodynamics, the emergent properties of the collective are considered. Just as 
it is possible to understand a thermal collection of atoms from scant features, such as the 
temperature or the entropy, so too can chaotic cosmologies be understood from a fractal 
dimension or a topological entropy. 

The singular cores of collapsing stars and the big bang are suspected to tend toward 
chaos Beyond conjecture^, attempts to conclusively identify chaos near singularities 

stirred debate |||6|]. The debate has centered on the mixmaster universe, an archetypal 
singularity. In the mixmaster model, the three spatial dimensions oscillate anisotropically 
out of the big bang and finally toward a big crunch @ . 

Relativistic chaos has the unique difficulty of demanding observer independent signa- 
tures. Many of the standard chaotic indicators such as the Lyapunov exponents and as- 
sociated entropies are observer dependent. The Lyapunov exponent quantifies how quickly 
predictability is lost as a system evolves. The metric entropy quantifies the creation of in- 
formation as time moves forward. They are both tied to the rate at which a given observer's 
clock ticks. The relativism of space and time rejects the notion of a prefered time direc- 
tion. In a curved space, the Lyapunov exponent and metric entropy become relative as well. 
Observer independent tools are needed to handle chaos in relativity. 

The challenge of relativistic chaos is well demonstrated in the debate over the chaoticity of 
the mixmaster universe. Barrow studied the chaotic properties of a discrete approximation 
to the full dynamics known as the Gauss map @||. Belinski, Khalatnikov and Lifshitz 
developed the Gauss map which evolves the universe toward the singularity from bounce 
to bounce of the axes ||. Barrow found the Lyapunov exponent, the metric entropy and 
topological entropy of the map. The dispute began when numerical experiments run in 



different coordinate systems found the Lyapunov exponents vanished [10-14|. Furthermore 



the Gauss map itself corresponds to a specific time slicing. The ambiguity of time was 
manifest . 

To conclude that the mixmaster approach to a singularity is indeed chaotic, an observer 
independent signature must be uncovered. As has been promoted elsewhere [[T5|-|ITf , fractals 
in phase space are observer independent chaotic signatures. Their existence and properties 
do not depend on the tick of a clock or the worldline of an observer. As well, fractals have 
an aesthetic appeal. The chaos is consequent of a lack of symmetries. The loss of symmetry 
in the dynamics is appeased by the emergent symmetry of the self-affine fractal [IS] in phase 
space. 

We exploit the observer independence of the fractal to show unambiguously that the 



Direct 2 and 3 dimensional numerical simulations of inhomogeneous cosmologies and gravitational 
collapse have so far failed to conclusively support or refute this conjecture Q. It is possible that 
numerical discretisation may suppress chaos as it causes a coarse graining of phase space similar 
to that found in quantum mechanics Q. 
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mixmaster universe is chaotic. This result was announced in Ref. Jl7| . 

Fractals in the minisuperspace phase space of the full dynamics are uncovered. Probing 
the full dynamics is always numerically intensive. To open a window on the numerical results 
we also study the Farey map, a discrete approximation to the dynamics related to the Gauss 
map. We focus on the map analytically but the real power of the conclusions is in the full 
dynamics where no approximation is made. Just as with a geographical map, the Farey 
map is used as a guide to navigate through the full dynamical problem. We use the map to 
locate the skeleton of the chaos and then verify its existence with the numerics. 

The bare bones of the chaos is revealed in an invariant subset of self-similar universes 
known collectively as a strange repellorQ. The strange repellor [F2U] is a fractal, nonattracting, 
invariant set (§|H). More familiar are fractal, attracting, invariant sets or strange attractors. 
Though the repellor is a tiny subset of all possible universes, it isolates the essential features 
of the system ||21|| . A typical universe will accumulate chaotic transient eras as it brushes 
past the strange repellor in phase space. 

The repellor is multifractal as shown in §11 A| . The Farey map divides the fractal into 



two complementary sets. A particularly elegant quality of the Farey map is the connection 
between these complementary sets and number theory. As elaborated in §11 B| , the multi- 



fractal can be completely understood in terms of continued fraction expansions and Farey 
trees. 

In addition to the usual multifractal dimensions, we stress the importance of the Lya- 
punov dimension ( §111 A]) . This dimension is built out of a coordinate invariant combination 
of the Lyapunov exponent (§ [III C|) and the metric entropy. Both the Lyapunov exponent and 



the metric entropy can then be reinstated as valid tools even in curved space. We introduce 
a method for handling the repellor as a Hamiltonian exit system, to facilitate calculation of 
the Lyapunov dimension ( §111 Bp . 



Maps have been put to good use before. We devote § [TV] to a detailed comparison of 
the Gauss map and the Farey map. There are two main differences we draw out. One 
distinction is in the character of the maps themselves. They are not topologically equivalent 
and therefore have different topological features such as topological entropy. The other 
notable difference is the technique for handling the maps. Previous methods treat the set 
of all universes. We concentrate on the repellor subset. The previous techniques and ours 
are complementary. 

Having isolated the repellor in the map, we use this insight as an x-ray to expose the 
chaotic skeleton in the full, unapproximated dynamics (§|V]). 

I. THE MIXMASTER DYNAMICS 

The richness of Einstein's theory is revealed in the spectrum of solutions to the dynamical 
equations. To solve these equations universe by universe, symmetries of the Hamiltonian are 



2 Page appears to describe a strange repellor in the minisuperspace phase space of a scalar field 



cosmology [19|. Remarkably, this paper was written in the same year strange repellors were first 



being described in chaos theory [20]. 
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sought. The entire universe must necessarily have a conserved Hamiltonian. Energy cannot 
leak into or out of spacetime. Since the Hamiltonian effects translations in time, symmetries 
of the Hamiltonian lead to invariants or integrals of the motion and hence solutions to the 
equations of motion. However, strong gravity is nonlinearly entangled and we cannot expect 
the universe to typically offer up symmetries of the motion. Einstein's equations are more 
likely than not to be fundamentally nonintegrable, and so chaotic. 

The singularity structure of the equations can indicate if the system is integrable. It has 



been shown that the mixmaster equations fail the Painleve test ||22|| . This suggests, but does 
not prove, that the system is nonintegrable. 

The fractals we find in the minisuperspace phase space show conclusively the emergence 
of chaos. To begin, we use dynamical systems theory on sets of universes and consider the 
minisuperspace of all possible mixmaster cosmologies. The minisuperspace Hamiltonian is 

H = (lna)'(ln&)' + (lna)'(hic)' + (ln&)'(lnc)' 

' a 4 + 6 4 + c 4 - 2(&V + aV + 6 2 a 2 )l (1.1) 



+ 4 



where a, b, c are the scale factors for the three spatial axes. For a full description of the 
geometry see Refs. |7|,|23|. The Hamiltonian constraint equation requires H = 0. The 



equations of motion governing the behaviour of the three spatial axes are 

(In a)" = ~ [(b 2 -c 2 ) 2 -a 4 ] et eye. (a, 6, c) . (1.2) 

A prime denotes d/dr where r is related to the cosmic time t through dt = abedr. The 
universe comes out of the big bang with two axes oscillating between expansion and collapse 
while the third grows monotonically. The overall volume of the universe expands to a 
maximum and then collapses. On approach to the big crunch, two spatial dimensions will 
oscillate in expansion and collapse while the third will decrease monotonically. Eventually 
a bounce occurs at which point the axes permute and interchange roles. 

The scale and expansion factors can be reparameterized in terms of the four variables 
(u, v, ||. As was done in reference ||, the initial conditions can be fixed on the surface 
of section a Q = 1, (hi.a)' > 0, without loss of generality 

in n v 3(1 + m V 

lna G = , ( ma )o 



'l + u + ul) 



, , 3fi G , -3u u 

\nb = — ■ -, (ln0) o 



[l + v Q + u v ) ° (1 + u + ul) 

3(v + u v )tt , 3(w + u 2 Vo 

lnCo= 7T^ , V' ln c = — ■ — z\. 1.3 

(1 + v + u a v ) {l + u + u 2 ) 

The relative sizes and velocities of the axes respectively are accounted for by v and u. The 
overall scale factors Q = (ln(afec))/3 and u — Q! grow monotonically away from the big bang 
to the maximum of expansion and then shrink monotonically toward the big crunch. 

The evolution of a given universe can be viewed as a particle scattering off a potential in 
minisuperspace. The right hand side of Eqns. (|1.2|) play the role of a potential. When this 
potential is negligible, the evolution mimics a simple Kasner model. The Kasner metric is 
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ds 2 = -dt 2 + t 2pa dx 2 + t 2pb dy 2 + t 2pc dz 2 



(1.4) 



The Kasner indices can be parameterized as 

1 + it 



Pa[U) 

Pb{u) 



— It 



{l + u + u 



2V 



u + u 2 



Pciu) = (iT^T^) ' (L5) 

and satisfy J2i=iPi = 1 = ELiPf- Moreover, p a (l/u) = p c (u), Pb(l/u) = Pb(u) and 
p c (l/u) = p a (u). As is customary, we take u to lie in the range [l,oo] to ensure a defi- 
nite ordering of the Kasner exponents: 

Pb<Pa<Pc- (1-6) 

In terms of conformal time, 

(In a)' ~ p a a oc exp(p a r) . (1.7) 

When the potential gains importance, the trajectory is scattered and enters another Kasner 
phase. To study the effect of scattering it is a good approximation to consider only the change 
on the Kasner indices of Eqn. ( |1.5| ). The full dynamical problem can therefore be reduced 
to a map which evolves the parameter u forward in discrete intervals of time. Belinski, 
Khalatnikov, and Lifshitz (BKL) reduced the dynamics to the 1-dimensional Gauss map. 
The Gauss map evolves u from bounce to bounce as the trajectory strikes the potential. We 
focus instead on a 2-dimensional map, the Farey map, which includes not only bounces but 
also oscillations Hj9| j2^j25[] . The Farey map has the additional property of being invertible 



and so preserves the time reversibility of Einstein's equations 



II. THE FAREY MAP AND THE STRANGE REPELLOR 

With the discrete time map we are able to isolate the repellor. The repellor is the 
chaotic subset of self-similar universes. It is the Hamiltonian analogue of the better known 
strange attractor. In dissipative systems, the volume of phase space shrinks and trajectories 
are drawn onto an attracting set as energy is lost. Since the universe is a self-contained 
Hamiltonian system, there can be no dissipation and so no attractors. Still, the notion of 
an invariant set can be an incisive characterization of the system. 

Nonattracting invariant sets are referred to as repellors since all trajectories which con- 
stitute the set are unstable in at least one eigendirection, and in that sense repel any near 
neighbours. The phase space volume is conserved as it squeezes in one direction while it 
expands in the other as it careens off the repellor. The transient chaos experienced by 
a typical, aperiodic universe reflects the passage of that trajectory near the core periodic 
orbits. 
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Khalatnikov and Lifshitz [|IJ were the first to introduce the u parameterisation and derive 
the evolution rules: 



u — > u — 1 if u > 1 
u — > 1 /it if u < 1 . 



(2.1) 



As it stands, this prescription is discontinuous at u = 1. A continuous map can derived by 
considering either of the two double transformations 



u 



u 



and 
and 



u 



u 



1/u' or 
u'-l. 



u-1 
= 1/u 

This leads to the two entirely equivalent continuous maps 

u n — 1, u > 2 
= < — , 1< U < 2 



(2.2) 



(2.3) 



and 



i/ 



n+l 



11,, 



■it., 



it,. 



li > 1 

< u < 1 



(2.4) 



Both choices have shown up in the literature over the years, but these days the first choice 
is prefered as the transitional step where u < 1 is hidden in the double transformation. 
By hiding this step, a definite Kasner ordering can be maintained. It is worth mentioning 
that the two maps use exactly the same u variable, even though the first can be formally be 
recovered from the second by making the substitution u = u — 1. 

By combining the w-map with its inverse in terms of v, one arrives at the two dimensional 
Farey map. We call the map a Farey map due to its connection to the Farey trees of number 
theory as discussed in § [11 A| . The Farey map evolves the two parameters forward by discrete 



intervals in time, F(u n , v n ) 



w n +i,f n +i)> according to the rule 
1, v + 1 , u > 1 (oscillations) 



u 



F(u,v) 



1 



1 

u — 1 V 



u < 1 (bounces) 



(2.5) 



The map describes the evolution of a universe through a series of Kasner epochs. If u > 2, 
then two axes oscillate in expansion and collapse while the third coasts. The axes will 
oscillate until u < 2 at which point a bounce occurs, the three axes interchange roles, and 
u is bounced about the numberline. In some studies, what we call oscillations are called 
epochs and what we call bounces are called eras. 

Running in parallel with the Farey map is the Volume map V(Q, u) described by 



n 



n+l 



a. 



1 + 



1 + u n + ui 



1 - 



u n {l + v n + u n v n ) 
2u r , 



(2.6) 
(2.7) 
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Like the Farey map, the Volume map is a good approximation to the exact evolution away 
from the maximum of expansion where uo flips sign. During the collapse phase r — > — oo, 
u) > and Q < 0. A nice feature of the combined Farey- Volume map, M(u,v,Q,u), is the 
way it splits the dynamics into oscillatory and monotonic pieces. 

Clearly, the monotonic behaviour of the l/-map excludes the possibility of periodic tra- 
jectories for the M-map. Since all chaotic behaviour relies on the existence of unstable 
periodic trajectories, this would seem to rule out chaos in the mixmaster dynamics. Indeed, 
similar reasoning has been used to claim that the full mixmaster equations cannot be chaotic 
|)2 7|| . The flaw in this reasoning was to neglect the non-compact nature of the mixmaster 
phase space. For example, the non-compact motion along a spiral of fixed radius appears as 
circular motion in an appropriate comoving frame of reference. Although a spiral trajectory 
never returns to its starting point, the motion is periodic in a projected subset of phase 
space. Similarly, the mixmaster dynamics allows unstable periodic trajectories in terms of 
the (u, v) coordinates. A projected strange repellor is all we require for chaos to occur. 
Consequently, we can drop the overall scaling of the V^-map and concentrate on the Farey 
map. Similarly, when we turn our attention to the full dynamics we will look for periodic 
behaviour in scale invariant quantities such as (In b) /(In c). In this way we effectively project 
out the regular evolution of the system. 

A universe initiated with some condition (u Q , v Q ) will evolve through a series of oscillations 
and bounces as prescribed by the map. For a typical universe, u and v will randomly be 
tossed about the numberline. An important subset of universes are formed by the orbits 
periodic in (u,v). The existence of a repelling, denumerably infinite, everywhere dense set 
of periodic orbits was first recognised by Bogoyavlenski and Novikov |28|| . Physically, the 



periodic orbits correspond to discretely self-similar universes. After one orbital period, the 
proportionality of the scale and expansion factors recur, while the volume of the universe has 
contracted overall. These collection of periodic orbits in (u, v) form a multifractal strange 
repellor. 

To isolate the strange repellor we search for a chaotic, nonattracting, invariant set. In 
this section we justify each term in this definition. 

Invariant set: The set of points which are invariant in time are located by the condition 

F k {u h Vi) = (ui,Vi) (2.8) 

for k = 1, . . . , oo. These fixed points lie on periodic orbits of period p < k. Iteration of the 
map fc-times, just returns a point on a period k orbit to itself. For Hamiltonian systems it 
is sufficient to consider the future invariant set as time reversal invariance can be used to 
find the complete strange repellor. For the Farey map the future invariant set is defined by 
the condition 

F k (ui) —Ui, (v arbitrary). (2.9) 

Consider for illustration the period 1 orbit. The invariant point satisfies F(u{) = u±. The 
map can be split into an oscillation map O and a bounce map B; 

0{u n ) = u n -l u n >2 (2.10) 

B{u n ) = — !— - u n < 2 . (2.11) 
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The equation F(ui) = u\ yields two possibilities: 0{ui) = u\ or B{ui) = u\. Only the 
bounce fixed point, which generates the equation 

«! = —!— , (2.12) 
U\ — 1 

has a solution. The solution is the golden mean 

% = ^ • (2.13) 

The occurrence of the golden mean is no accident. As we discuss in the next section, it 
is the first leaf on the irrational Farey tree. No matter how far into the future we try to 
evolve this initial condition, we get u% back. A universe initiated with axes and velocities 
parameterized by u\ will scale in a self-similar manner as u is returned to u\ with each jump 
forward in time. 

Another particularly simple set of fixed points can be located. These are the maximal 
value of u along a period k orbit. The maximal value of u corresponds to the largest number 
of oscillations before a bounce. The equation O k ~ l B{u) = u generates 

u k (u k -k) = l, (2.14) 

whose solution is 

u k = • (2.15) 

These are the silver means (and can be found among the kth order leaves in the Farey tree). 
As elaborated in the next section, all of the fixed points comprise a countably infinite set 
of irrationals with periodic continued fraction expansions described by Farey trees. In this 
sense, we know all the fixed points. 

Nonattracting: We can verify that all of the periodic orbits are unstable and so are not 
attracting but rather repelling. The repellor is the intersection in phase space of the unstable 
and stable manifolds. Since it is unstable in one direction, a near neighbour to a periodic 
orbit will deviate off the invariant set. Conservation of the Hamiltonian requires the other 
eigendirection be attracting. 

In the 2-D phase space defined by (u, v), we show here that the unstable manifold is 
composed of fixed points along u. This is the future invariant set given by Eqn. ( |2.9|) . Since 
time reversal corresponds to (u, v) — > (v,u), the stable manifold is made up of the fixed 
points along v. The stable and unstable manifolds intersect at each point along a periodic 
orbit. The collection of such points form the strange repellor. Since they are invariant, any 
point which is on both the unstable and stable manifolds will be mapped to points which 
by necessity inhabit at the intersection of stable and unstable manifolds. 

To demonstrate explicitly that the periodic orbits are unstable in the u direction but 
attracting in v consider an initial trajectory in the vicinity of a period p orbit. We look 
at the evolution in the u direction first. The aperiodic universe begins with u Q = u Q + 5 
where u a is the fixed point along the periodic universe. A period later, the deviation from 
the periodic orbit has grown to 8 P defined by 
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F p {u ) = u + 5 P = F p {u + S ) . (2.16) 
Expanding the right hand side we can relate the evolved deviation to the initial deviation 

S p = CpS (2.17) 

where the stability coefficient is 

dF p {u n 



P dUr. 



= F'{u )F\u l )...F\up^ 1 ) . (2.18) 

The elements in the final product are evaluated at the points along the unperturbed orbit. 
The periodic orbits are stable if the magnitude of S p shrinks and unstable if the magnitude 
grows. In terms of the stability coefficient the periodic orbits are stable if \c p \ < 1 and 
unstable if \c p \ > 1. Taking the derivative of the map 

1, 1 {u > 2) 

\F(u,v)'\ = \ l l . (2.19) 

1 7^2 > ~ (u<2) 

Along any periodic orbit, at least one point must fall between 1 and 2. Let u m \ n be the 
smallest value of u along the orbit. It then follows from ( |2. 18|) that 



\c P \ > (u min - I)' 2 > 1 . (2.20) 

Therefore all orbits on the repellor are unstable in the u eigendirection. Any universe 
which begins near an exact periodic orbit eventually deviates away from that orbit in u 
and is not attracted onto it. Since a typical universe is confined to bounce around the 
numberline forever it will eventually scathe past a periodic orbit before being repelled off 
only to accidentally stumble onto the repellor again. Thus, a typical universe will scatter 
around intermittently hitting chaotic episodes as it jumps on and off the repellor. 

Repeating the stability analysis in v , we find \c\\ = 1/v 2 < 1. Near neighbours tend to 
follow periodic orbits in the v direction. The repellor along u is an attractor along v. 

Chaotic: We have so far established that the periodic orbits comprise a nonattracting, 
invariant set. We can demonstrate that this set is chaotic by showing it has a positive 
topological entropy. The topological entropy, in analogy with the thermodynamic entropy, 
measures the number of accessible states on the repellor. The number of states on the 
repellor is equivalent to the number of fixed points so 

H T = lim y\nN(k) (2.21) 

where N(k) is the number of fixed points. For a nonchaotic set, the number of fixed points 
is either finite or grows as a finite power of k, so = 0. 

The fixed points at order k are found by solving equation (|2]9|): F k {ui) = Ui. To count 
the number of fixed points we can count the number of such possible equations. For a given 
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period p < k, F k can be broken into O m oscillations and B k ~ m bounces. Since O and B do 
not commute, the order in which they occur leads to different possible solutions for the ttj. 
The number of ways to combine m O's and k — m B J s is 

k \ rn\ . 

2.22 



m J m\{k — m) 

For example, after k = 4 iterations of the map and m = 2 oscillations and k — m = 2 bounces 
there are 6 possible permutations: 

OOBB{ui) = Ui 
BBOO{ui) = Ui 
BOOB{ui) = Ui 
OBBO(ui) = v,i 
OBOB{ui) = Ui 

BOBO(ui) = Ui (2.23) 

However, the first 4 are all cyclic permutations of OOBB. Cyclic permutations must lie 
along the same orbit. Therefore, the solutions to the first 4 equations yield the points along 
the same period 4 orbit. Similarly, the last 2 are cyclic permutations of each other. They 
represent the recurrence of the two points along the period 2 orbit. 

The total number of points belonging to the future invariant set with period p < k is 
given by summing over all possible combinations of oscillations and bounces: 

- k 



m=0 



There are thus 2 k — 1 words of length k that can be built out of a 2 letter alphabet. The 
topological entropy is then 

#£ = ln2 . (2.25) 

This entropy is independent of phase space coordinates. The topological entropy for the 
full 2-D map F(u,v) will be twice this quantity as the strange repellor is formed by the 
intersection of 2 k — 1 horizontal and vertical lines. Thus, there are (2 k — l) 2 roots of Eqn. ( |2.8| ) 
at order k and Ht = 2 In 2. The result = In 2 was first found by Rugh |10| using symbolic 
dynamics to describe typical, aperiodic orbits in u. 

The topological entropy will be the same for all maps which are topologically conjugate; 
that is, which can be related by a continuous, invertible, but not necessarily differentiable, 
transformation of coordinates. In particular, this tells us that the map F(u) can be obtained 



from the shift map, horseshoe map or generalised baker's map [29fl. We mention that the 
closely related Gauss map has a much higher topological entropy = 7r 2 /(61n2). As 
discussed in §[TV], the two maps are not topologically conjugate. We do not expect then for 
their entropies to be the same. 

The bare bones of the chaotic scattering has been exposed in the repellor. Taking this 
skeleton, we can now show that the repellor is in fact strange, that is, fractal. 
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A. The Repellor is a Multifractal 



We can create an illustrative picture of the fractal set of self-similar universes. Distribut- 
ing the fixed points u along the numberline, the collection of periodic points form a fractal in 
phase space. A fractal is a nowhere differentiable, self-affme structure. It cannot be undone 
by a coordinate transformation. The emergence of a fractal distribution of invariant points 
is an observer independent declaration of chaos. 

For multifractals there is information in the way points are distributed. There is an 
underlying fractal structure with an architecture of distributed points built on top. The 
fractal in phase space is simply constructed by locating all of the points along period p < k 
orbits numerically and plotting them in a histogram as is done in Fig. 1. The histogram 
was generated by solving for all roots of Eqn. (|2~9l ) up to and including k = 16. The 
histogram reveals how points are distributed in the future invariant set, otherwise known as 
the unstable manifold of the strange repellor. The self-similarity of the distribution is clear. 

The distribution of points can be understood in terms of the period of the orbit. As k 
increases, the lowest period orbit, namely the p = 1 orbit, is visited the maximal number 
of times, that is, k times. On the other hand, the maximum value of u at order k lies on a 
p = k orbit, and so is visited only once. 

Points on the repellor are clustered in the interval 1 < u < 2. The combinatorics of 
building words out of a 2-letter alphabet (O, B) of oscillations and bounces favours the 
small u tower. Consider the integer interval [n, n + 1] . A root in this interval corresponds to 
the sequence, or word, O^-^B occuring somewhere along the orbit. The complete orbit is 
a sentence of O l B words, eg. (0 2 B)(B)(0 5 B)(OB), repeated in a cyclic fashion. Since the 
number of n-letter words that can be formed from a 2-letter alphabet is 2 n , it follows that 
the fraction of roots in each interval is 

Pn = 2~ n (2.26) 
Note that the distribution is correctly normalised since 

oo 

£pn = l. (2.27) 

n=l 

The exponential fall off in the density of points on the repellor is clearly evident in Fig. 1. 

The parameterization of the axes in terms of u in Eqn. ( pL.5| ) shows that the small u 
orbits typical of the repellor are those with axes of similar scale and speed. These universes 
have axes which frequently switch from expansion to collapse. It follows that the strange 
repellor corresponds to the most isotropic mixmaster trajectories possible. 

In the full 2-D phase space of the Farey map, F(u, v), the unstable manifold appears as 
a forest of vertical lines, while the stable manifold appears as a forest of horizontal lines. 
A portion of the future invariant set, the unstable manifold of the repellor, is displayed in 
Fig. 2. The collection of lines appear to form a multifractal cantor set. To confirm this, we 
calculate the fractal dimension of the set. 

The fractal dimension can be thought of as a critical exponent. The dimension is a 
measure of the length of the collection of points. To find its box-counting dimension, cover 
the set with N(e) boxes of size e. As the size of the boxes is taken infinitesimally small, the 
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number of boxes needed to cover the set grows. For fractals, the number of boxes needed 
grows faster than the scale shrinks. The structure becomes more and more complex as 
smaller and smaller scales come into focus. The critical exponent can be defined as the 
number which keeps 

lime Do iV(e) (2.28) 

finite. Any exponent greater than D Q would result in zero length, and an exponent smaller 
would result in infinite length. In this sense D a is a critical exponent. 

This dimension can be generalized to include not just the scaling property of the fractal 
but also the distribution of points on top of the underlying foundation. This leads to a 
continuous spectrum of critical exponents or dimensions. The spectrum of dimensions is 
more commonly expressed as 

fl , = ^ h !!ggM!, (2 . 29) 

q q - 1 e-a In e v ' 

where N(e) are the number of hypercubes of side length e needed to cover the fractal and 
Pi is the weight assigned to the i th hypercube. The p,'s satisfy Ya=x Pi — 1- The standard 
capacity dimension D Q is recovered when q = 0, the information dimension when q = 1, the 
correlation dimension when q = 2 etc. For homogeneous fractals all the various dimensions 
yield the same result. The multifractal dimensions D q are invariant under diffeomorphisms 
for all q, and D\ is additionally invariant under coordinate transformations that are non- 



invertible at a finite number of points [ 30 1 



The quantitative importance of the fractal dimension is expressed in terms of the final 
state sensitivity f{5) |3l|. This quantity describes how the unavoidable uncertainty in 



specifying initial conditions gets amplified in chaotic systems, leading to a large final state 
uncertainty. By identifying a number of possible outcomes for a dynamical system, the 
space of initial conditions can be divided into regions corresponding to their outcome. If 
the system is chaotic, the boundaries between these outcome basins will be fractal. Points 
belonging to the basin boundary are none other than the chaotic future invariant set. The 
function f(5) is the fraction of phase space volume which has an uncertain outcome due to 



the initial conditions being uncertain within a hypersphere of radius 6. It can be shown pi 
that 



f(S) ~6 a , a = D-D , (2.30) 

where D is the phase space dimension and Dq is the capacity dimension of the basin bound- 
ary. For non-chaotic systems a = 1 and there is no amplification of initial uncertainties, 
while for chaotic systems < a < 1 and marked final state sensitivity can occur. 

It can be argued that the capacity dimension must equal the phase space dimension for 
the mixmaster repellor. The periodic orbits are given by the periodic irrationals as explained 
in §|IB[ The periodic irrationals are dense on the numberline. Formally this means there is 
always another periodic irrational e away from a neighbour for all e > 0. Since the periodic 
orbits are dense, there will always be an infinite number of fixed points in any box. It is 
sufficient to cover the set then with N(e) ~ 1/e 2 boxes. Therefore, the basic box counting 
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dimension is D a = D = 2 |[32|| , though it converges slowly. In the infinite time limit this can 
be interpreted as an ultimate loss in predictability since f(5) = 1. The mixmaster is very 
mixed. 

While the box counting dimension saturates at the phase space dimension, the more 
heavily weighted dimensions of Eqn. fl2.29|) do not. Taking the information dimension of the 
fractal in Fig. 2, we find 

D\ = 1.87 ±0.01 . (2.31) 

Fig. 3 shows the fitf\ used to determine D\. The quality of the fit gives us confidence that 
the strange repellor is truly fractal. 

By forming the intersection of the Farey map's stable and unstable manifolds, ie. by 
solving Eqn. Q2.8|) for all fixed points (u,v), we uncover the strange repellor. A portion of 
the strange repellor is shown in Fig. 4 using all roots up to k = 12. While it is possible to 
find the information dimension of the strange repellor directly from Fig. 4, the effort can 
be spared. Since the Farey map is Hamiltonian, we know the stable and unstable manifolds 
of the strange repellor share the same fractal dimensions. That is, D\ = D\ = 1.87 ± 0.01. 
Now, the repellor is formed from the intersection of the stable and unstable manifolds, so 
its information dimension is simply 

D 1 = D\ + D V X - D = 1.74 ±0.02. (2.32) 

Since D Q ^ D\, we have confirmed the multifractal nature of the strange repellor. The 
dimension of Eqn. ( |2.32j ) is calculated with all fixed points occurring at order k — 16. As 
k increases this number will continue to grow slowly as the rarified regions of the fractal 
continue to be populated. Thus, (|2.32|) actually represents a lower bound. However, we 
expect the ultimate value should only differ in the second decimal place. 



B. Farey trees 

The fractal of Fig. 4 reveals two complementary sets. The sequence of gaps correspond to 
the rational numbers. In a complementary fashion, the strange repellor is made up from the 
periodic irrationals. Both sets are of Lebesgue measure zero. This division of the numberline 
can be understood in terms of the properties of the continued fraction expansions (cfe). We 
can write any number, rational or irrational, in terms of a cfe. While the connection to cfe's 
was known and explored by Belinskii, Khalatnikov and Lifshitz, they explicitly ignored the 
periodic irrationals [EJ. We focus on the periodic irrationals as they constitute the strange 
repellor. 

Consider some initial condition, u a , not necessarily on a periodic orbit. We can decom- 
pose any number into its integer part and some left over: 



3 The fit is actually for the fractal formed by taking a horizontal cross section through Fig. 2. The 
final answer is obtained by adding 1 to this number. 
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u = m + x D , (2.33) 

where m = [u Q ] denotes the integer part of u Q and x Q the fractional excess. According to the 
map, the next value of u following a bounce is u\ = l/x Q . Now decompose l/x Q = mi + xi, 
so that mi > 1. Solving this for x Q we find x Q = l/(mi + £i). At the next bounce «2 = 
so that Xi = l/(m2 + £2)- In this way we generate the cfe for u Q 

u = m A ; (2.34) 

mi H y 

m 2 H — 

where m^ > 1 Vi. The integers — 1 represent the number of oscillations between bounces. 
In shorthand form the cfe can be written as u Q = [m , mi, m 2 , m 3 , m 4 , . . .]. 

The map naturally distinguishes numbers on the basis of their cfe. A rational number 
can be written as the ratio of integers x Q = p/q. Consequently, the cfe is finite. The fact 
that its finite means at some iteration u n+ i = 00 and rationals are tossed out of the map. At 
the opposite extreme, the periodic orbits on the repellor have infinite cfe's which repeat. For 
example, an orbit of the form {O l ^ 1 BO^ l B} 1 where the curly brackets denote a repeated 
pattern, has the expansion 

u = i + -j . (2.35) 

* + 1 

% + — 



1 

3 + 



i + ... 

In shorthand form this reads u Q = [{i,j}]. Lagrange showed that the necessary and sufficient 
condition for an irrational number to have a periodic cfe is for it to be the root of a quadratic 
equation with integer coefficients. It is easy to prove that the Farey-map gives such equations 
for the periodic orbits at every order in k. 

The map generates a Farey tree Consider the initial condition u Q = 1. With the first 
iteration of the map, ie. with the first attempt to evolve that universe forward in time, u is 
thrown to infinity. After two iterations, a universe with the initial condition u Q = 2 will be 
thrown to infinity. At the kth application of the map, after the universe has evolved forward 
in k jumps, another group of universes whose initial conditions were rational numbers are 
discarded. The pattern of discarded terms is a Farey tree |3"3fl. Farey trees also arise in the 



quasiperiodic route to chaos described by the circle map |33 
In the interval [1, 2] the Farey tree can be written as 




(2.36) 



The bracketed terms each comprise a level of the Farey tree. For all k > 2, there are 2 h ~ 2 
leaves at order k. Every rational number in the interval [1,2] occurs exactly once somewhere 
on the Farey tree. Each leaf on the Farey tree has a continued fraction expansion that 
satisfies 
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J2rrii = k. (2.37) 



i=l 



Now, we see that each m, corresponds to a m^-letter word 0^ mi ~^B. Each rational number 
corresponds to a sentence of these words. For example, at level k = 5 we have the Farey 
number 10/7 = [1,2,3], which gives rise to the sentence BOBOOB . At the same order we 
also have 6/5 = [1,5] = BOOOOB, and 9/5 = [1,1,4] = BBOOOB etc. In other words, 
a universe which began with u = 10/7 will follow the pattern of BOBOOB and then gets 
tossed to infinity. The universe then evolves with two axes fixed and the third expanding as 
oc t. Such a universe has the form of a Rindler wedge xR 2 |35[| . 

At the opposite extreme from the rationals which escape the map are the periodic subset 
of irrationals. This explains the occurrence of the golden mean as the period 1 orbit with 
no oscillations. Its cfe is 



[{1}] . (2.38) 



1 



1 + ... 



Numbers with a periodic cfe are a set of measure zero among the irrationals, though they 
are dense on the numberline. 

An irrational Farey tree for the repellor can be constructed using the rational Farey tree as 
a seed. Consider the construction in the interval [1,2]. For each leaf on the rational Farey tree 
given by [1, mi, m^, ... , mt~), add a leaf to the irrational tree given by [1, {mi, m 2 , . . . , m^}]. 
There are 2 times as many leaves on the irrational Farey tree than there are on the rational 
Farey tree. This follows since all rationals have two continued fraction expansions. The two 
cfe's differ in how they end. In general 

[m , ...,m j + l] = [m , . . . , rrij, 1] . (2.39) 

Consider 10/7 = [1, 2, 3]. We see that 10/7 also equals [1, 2, 2, 1]. Thus, 10/7 « 1.429 has two 
corresponding leaves on the irrational Farey tree, (Vl~5 — l)/2 « 1.436 and (V85 + 5)/10 ~ 
1.422. Consequently each rational generates two nearby irrationals. 
A number theory result tells us that 



V 

a± 

q 



(2 - 40) 



where a± stands for the two irrational numbers that lie on each side of their rational seed. 
As k — > oo the leaves on the rational and irrational Farey trees get closer and closer together. 
Another result from number theory tells us that a + > p/q and a_ < p/q (hence the names). 
The 1/q 2 gaps around the rational Farey numbers seen in Figs. 1&2 are a consequence of 
this number theory result. 

The asymptotic distribution of roots goes as follows: At level k there are 2 k ~ 2 leaves on 
the rational Farey tree in the interval [1,2]. Each leaf produces 2 leaves on the irrational 
Farey tree. Thus, there are 2 k ~ 1 roots in the interval [1, 2] at level k. This matches our 
earlier results that there are ~ 2 k roots at order k, half of which lie in the interval [1,2]. 

Our brief excursion into number theory has produced a neat picture: The rational Farey 
tree gives us all the numbers that get mapped to u = oo, while the irrational Farey tree 
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gives us all the periodic orbits. The leftovers are all the irrationals, save the rational and 
irrational Farey trees, both of which have measure zero. Thus, all typical trajectories are 
aperiodic, unbounded and of infinite length. These are the trajectories typically studied in 
the literature. Our study is complementary. 



C. Lyapunov Exponents 

In flat space, the extreme sensitivity of the dynamics can be quantified by Lyapunov 
exponents and the related metric entropy. The Lyapunov exponents determine how quickly 
in time trajectories diverge. The metric entropy measures the rate at which information 
is created. Since the exponent and the entropy are rates, they connect directly with the 
rate at which time pushes forward. Clearly, the observer dependence of the rate at which 
clocks tick make these tools suspect in the hands of relativists. Relativists abandon such 
coordinate dependent indicators. 

Although Lyapunov exponents are observer dependent and therefore ambiguous in gen- 
eral relativity, we see in §[111 A| that they are related to an observer independent quantity; 



namely, the Lyapunov dimension. Since there is still utility in them, we take the time to 
compute some Lyapunov exponents. 

Since we know all the periodic orbits of the Farey map analytically, we are in a position 
to calculate the Lyapunov exponents for trajectories belonging to the Farey repellor. As a 
warm up we calculate the Lyapunov exponents for the the golden and silver mean orbits. 
Then we use the irrational Farey tree to write an analytic expression for the Lyapunov 
exponent for any periodic orbit on the repellor. 

In close analogy to the manner in which the stability coefficient was found, the Lyapunov 
exponent for a given orbit is defined by 

A = & iSM^MI (2-4i) 

n=0 

where the u n are the points along the orbit under scrutiny. The Lyapunov exponents along 
v are the negative of those along u as we now show. The points along a periodic orbit are 
the same whether time runs forward or backward. Time reversal corresponds to inverting 
the map. Note that (F~ r )'(u) = 1/F'(u) and therefore ln|F _1/ (M)| = — ln|F'(-u)|. Since 
F^iu^v) = F{y,u) it follows that ln|F'(u)| = —Yn.\F'(u)\. This identity in eqn (p. 41 



shows that 

\ v = -X u , (2.42) 

as must be the case to conserve the phase space volume in the Hamiltonian system. 
In illustration consider the golden mean period 1 orbit. All the u n = U\ and 

F\u n ) = -{u x - I)' 2 = -(mO 2 . (2.43) 

From Eqn. ( |2.41| ) it follows that 

Ai= lim I(T-l)lnn?, (2.44) 

1 — »oo ± 
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which reduces to 



Ai = 2 In 1 ~ 0.9624 . (2.45) 

Similarly, for the silver means ( 2.15|) we can find the exponent. We know the maximum 
value of u along the silver orbit is given by w max = [\Jk 2 + 4 + k)/2, so the sum of the ln's 
in Eqn. ( |2.41| ) becomes the In of the product, 

F'(u )F'(u 1 )...F'(u k ^ 1 ) = -(u max )- 2 . (2.46) 

After k applications of the map, F k (u) = u, so that after T application, a given value of u n 
has repeated [T/k] = T/k times. Thus T must be a multiple of k. The sum in the ( p. 41 
can be written 



T-i \T-l 
£ln|F'( Un )| = 

n=0 



21nw max (2.47) 



Which gives 



2, k+VW+4: . A , 

Xk= k ln 2 • (2 ' 48) 



We have found the Lyapunov exponents for the two extreme cases of the period 1 orbit 
and the longest period k orbit. We now use the irrational Farey tree to write down the 
general expression for Lyapunov exponents: 

A = 2(j>,) $>([{m i; m 2 ,...}]) , (2.49) 

\i=0 / eye 

where the second sum is taken over all terms in the cycle {mi, mz, ■ ■ ■}■ Since rrii > 1, each 
term in the sum of logs is greater than zero. The Farey tree has given a nice compact form 
with which to express all the Lyapunov exponents. From ( |2.49| ) it is easy to recover the 
results for the golden and silver means: 

2 rri ,. 2, k + ^/k 2 + 4 
X k = - \n[{k}} = - In V - . (2.50) 

We can also calculate the average Lyapunov exponent for the periodic orbits. Consider 
a very long orbit of length / ^> 1. On average, 1/2 points along the orbit will contribute 
to the sum in ( |2.49[ ) since half of the points will lie in the interval [1, 2]. Thus, the average 
Lyapunov exponent is given by 

A P =yE ln ([{"Hi m 2 , -..}]) =< ln([...]) > , (2.51) 

eye 

where <> denotes the average of all contributions. Since we know the probability that each 
rrii equals a given integer n goes as 2~ n , it follows that 
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A = nm Emi Em 2 •■■ Em, M[{ m i; ma, m ; }] ^ ^ 



Emi Em2 "' Em; Pi 

where the density pi is given by 

Pi = II Prm = 2 -( mi+m2+ - +m; ) . (2.53) 

i=l..l 

The excellent convergence properties of the sums in ( |2.52| ) ensure that a finite truncation is 
able to provide a good estimate. Using I = 6 and summing mi up to 30 and the 1712.. rriQ 
up to 10 we find A p = 0.793. We can check this result by numerically evolving a large 
number of periodic orbits. Using Maple to find all 2 15 — 1 = 32767 periodic orbits at order 
k = 15 and then evolving these orbits numerically, we find an average Lyapunov exponent 
of A p = 0.792 ± 0.002. This is in excellent agreement with the finite truncation of the exact 
sum. 

In contrast to the periodic orbits, typical aperiodic orbits have vanishing Lyapunov ex- 
ponents. The average Lyapunov exponent for a typical aperiodic orbit can be approximated 
by 

2 

<A>= — lnw av . (2.54) 

From Fig. 5 we see that w av — > 00 as the number of iterations of the map grows large. Thus, as 
first noted by Berger |13||, the average Lyapunov exponent for aperiodic trajectories tends to 



zero as n — > 00. This behaviour is characteristic of a chaotic scattering systems. Trajectories 
on the strange repellor have positive Lyapunov exponents while typical scattered orbits have 
Lyapunov exponents that tend to zero. The chaos in these systems is called transient as the 
brief chaotic encounter with the strange repellor is followed by regular asymptotic motion. 
Of course, all of these statements should be made with extreme care in general relativity as 
Lyapunov exponents are not gauge invariant. 



III. THE LYAPUNOV DIMENSION, TRANSIT TIMES AND AVERAGE 

EXPONENTS 

A. The Lyapunov Dimension 

We have extolled the virtues of fractals in phase space as coordinate independent signals 
of chaos. In this section we relate the fractal dimension of the invariant set to important 
dynamical quantities such as Lyapunov exponents and metric entropies. The Lyapunov 
dimension Dl combines these coordinate dependent quantities into an invariant combination. 

Remarkably, it has been shown that Dl equals the information dimension D\ for typical 
chaotic invariant sets. The relation D\ = Dl has been rigourously established for certain 
dynamical systems |36[] and has been numerically confirmed for many typical systems [57 



Specifically, it has been conjectured that p8|,|39| 



Aa+i 
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Here < r > is the lifetime of typical chaotic transients, A is largest integer such that Ai + 
A 2 .. + Aa > 0, and the Lyapunov exponents are ordered so that Aj > Aj+i. The equality of 
the information and Lyapunov dimension is believed to hold for most continuous dynamical 
systems, although a rigourous proof has only been given for discrete maps. 



A heuristic derivation of this result can be found in the text of Ref. ||29|| . We sketch that 
reasoning here for the 2-D Farey map. The repellor marks the intersection of the stable and 
the unstable manifolds. While it is repelling (unstable) in the direction of u, it is attracting 
(stable) in the direction of v. We can define a natural measure on the stable manifold, 
unstable manifold and on the set itself. Ordinarily, an initial condition will produce an orbit 
which leaves the repellor never to return. Consequently, the number of trajectories near the 
repellor decays with time. The measure on a set is loosely related to the number of points 
which hang around the set. Consider some number N(0) of random values for u scattered 
about the numberline. These orbits are evolved by the map and eventually expelled. After 
a large number of iterations the only trajectories that remain belong to the invariant set, or 
are at least very close to a trajectories belonging to the invariant set. Thus, the measure of 
the repelling set can be defined as 

fi(t) ~ N(n)/N{0) ~ exp(-n/ <r>) . (3.2) 

The time scale <r> characterizes the decay time of typical trajectories leaving the repellor. 
In other words, < r > is the lifetime of typical chaotic transients. The measure can also be 
connected to the notion of the length of the set and consequently to the fractal dimension 
through 

e 2 - D "~/i, (3.3) 
where e ~ exp(— <A> n) and the Lyapunov exponent is defined by 

I 

<A>= lim — — Y A^ (3.4) 



The information dimension of the unstable manifold, D™, makes an appearance in ( |3.3|) 
since the measure accounts for the distribution of points as well as their location along the 
numberline. Taking the natural logarithm of Eqn. (|3.3j ) yields 



D« = 2-— i . (3.5) 

<A><r> 

Thus, we can relate the information dimension of the unstable manifold (the periodic ir- 
rationals in u) to the positive Lyapunov exponent and decay time of the Farey exit map. 
The dimension of the repellor is the sum of the dimensions of the stable and the unstable 
manifolds. For the conservative and hence invertible Farey exit map, the dimensions of the 
two manifolds are equal. Thus, the information dimension of the strange repellor is simply 
D\ = 2D\ — D. Consequently, the Lyapunov dimension is 

Di 2 ( I ' 



<T>< A> 



2 ^L, (3,) 
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where h(/i) is the the metric entropy 



%)=<A>-— — . (3.7) 

<T> 

We see that Dl is given by the ratio of the metric entropy to the average Lyapunov exponent. 
Although neither Lyapunov exponents nor metric entropies are coordinate invariant, their 
ratio is. If a system is chaotic, that is, if D\ ^ 0, a coordinate system can always be found in 
which h and <A> are both finite and non-zero. From a dynamical systems perspective, such 
coordinate systems are preferable as they allow us to reinstate both Lyapunov exponents 
and metric entropies as useful chaotic measures. 

In order to implement Dl, we derive <r> and <A> for the mixmaster model in the 
following sections. 



B. Hamiltonian Exit Systems 

In our case, the system does not create an efficient repellor. Typical universes after being 
scattered off one periodic universe will eventually happen across another. The battle to toss 
trajectories off as they continually wash back ashore is constant. Because points which are 
discarded from the repelling set eventually return, the time <r> needed to discard typical 
trajectories from the repellor is infinitely long. 

In principle, there is nothing wrong with a repellor that is revisited by previously scat- 
tered trajectories. In practice, it is easier to handle systems where scattered trajectories are 
discarded once and for all. We introduce a method for turning our thwarted repellor into a 
cleaner Hamiltonian exit system in this section. 

The mixmaster dynamics presents some unique challenges as a dynamical system. It has 
some characteristics of a chaotic billiard but this picture is upset by the non-compact 
nature of its phase space. The mixmaster also shares characteristics of a chaotic scattering 
system, but this picture is upset by the lack of absolute outcomes. Only a very special set 
of initial conditions lead to universes which terminate after a finite number of bounces. 

There are other dynamical systems that combine features of chaotic billiards and chaotic 



scattering. These are known as Hamiltonian exit systems [fil |. A standard example is a 
chaotic pool table. Trajectories can bounce chaotically around the table before falling into 
a particular pocket. When the system is chaotic, the pocket a ball finishes up in depends 
sensitively on initial conditions. Each pocket has a basin of attraction in phase space, and 
the borders between the basins of attraction can be fractal. As discussed in § [11 A| , the fractal 
dimension of the basin boundaries provides a direct measure of the sensitive dependence on 
initial conditions. 

In many ways, adding exits to a chaotic billiard is the best way to study the dynamics of 
the original closed system. By opening exits we can expose the underlying chaotic invariant 
set that encodes the chaotic behaviour. The process would be familiar to an archeologist 
looking for bones. By running water through a sieve, dirt is washed away to expose the 
bones. Similarly, opening holes in a chaotic billiard allows most trajectories to escape 
leaving behind the skeleton of unstable periodic orbits. As in archeology, some care has to 
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be taken with the choice of sieve as bad choices can lead to the loss of bones along with the 
dirt. 

For the mixmaster system, the placement of the pockets follows naturally. There are 
already three infinitesimally thin pockets. The three trajectories which lead out to these 
three pockets correspond to the Rindler universes with Kasner exponents {p a = pb = 0, p c = 
1}> {Pa = Pc = 0, pb = 1} or {pb = p c = 0, p a = 1}. All we have to do is widen the pockets a 
little and the strange repellor will lie exposed. Since the three pockets correspond to large 
values of u, we know from our study of the Farey map that typical trajectories spend most 
of their time near a pocket. In addition, trajectories on the strange repellor are concentrated 
at small u values, and are thus far from the natural pockets. Only a small percentage of the 
strange repellor will be lost when we widen the pockets. 

A nice visual representation of the mixmaster billiard is provided by the effective potential 
picture. The minisuperspace potential is given by 

V = (a 4 + 6 4 + c 4 - 2a 2 b 2 - 2b 2 c 2 - 2c 2 a 2 )/(abc) . (3.8) 

The potential is best viewed by making the change of coordinates (a, b, c) — > (fl, 

Q = - ln(afec) , (3.9) 
3 

/3 + = il,g), (3.10) 

so that 

V = e- 8f3+ + 2e 4/3+ (cosh (4^3/3-) - 1) (3.12) 
- 4e" 2/3+ cosh(2v / 3/?-) . 

Equipotentials of V are shown in Fig. 6. The natural pockets occur at the accumulation 
points Q\ = 0, 62 = 27r/3 and #3 = 4ir/3. The angle 9 is measured from the (3 + axis. 
Mixmaster trajectories can be thought of as a ball moving with unit velocity as it bounces 
between the triangular walls. The walls are also moving, but at half the velocity of the ball. 
Equipotentials of V correspond to the position of the wall at different stages during the 
mixmaster collapse. 

Typical trajectories oscillate around in one of the corners for a long time before bouncing 
out to a new corner where they oscillate around for a long time and so on ad infinitum. These 
are the typical, aperiodic mixmaster trajectories. Thus, the mixmaster dynamics leads to 
three accumulation points at the corners of the triangle. This behaviour was first noted by 



Misner |7j , and was later studied numerically by Creighton and Hobill |42| . 

In addition to the aperiodic trajectories, there are two special classes of mixmaster trajec- 
tories. One class corresponds to the rational numbers. After a finite number of oscillations 
and bounces they take the perfect bounce and head straight out one of the infinitesimally 
thin pockets. The other class of special trajectories form the strange repellor and correspond 
to the periodic irrationals. These trajectories regularly bounce from corner to corner and 
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spend little time oscillating in each corner. Consequently, trajectories on the strange repellor 
are unlikely to exit the system when we widen the pockets. 

By widening the natural pockets, we create exits. For example, the first pocket becomes 
the angular region [-A9, AO]. As the collapse proceeds, Q — > — oo, and we can relate AO to 
the map parameter u via 



A6 = 2 arctan ( \f?> 
u \ 



u 



1 



u — 1 
1 1 

2u Au 3 



2tt 



(3.13) 



We see that choosing an exit value for u sets the angular width of the pockets. In the limit 
Wexit - > °o we recover the original mixmaster pockets with zero angular width. In Fig. 7 
we display equipotentials of the minisuperspace potential with pockets corresponding to the 
exit value M exit = 16. While the pockets appear quite wide, we see from Eqn. ( 2.26|) that they 



have little effect on the strange repellor. Few bones will be lost from the chaotic skeleton. 
For Mgxit = 16 less than 0.003% of the strange repellor will be lost out the corner pockets. 
The wider we make the pockets the faster the strange repellor is exposed. For numerical 
studies of the full dynamics we chose w ex it = 8, giving exits ~ 1.7 times wider than those 
shown in Fig. 7, but still small enough to ensure that less than 1% of the chaotic skeleton 
is lost. 

With the exits in place the mixmaster behaves like a classic chaotic scattering system. 
If we start Nq randomly chosen trajectories, there will be 



N(n) = N exp 



n 



<T> 



(3.14) 



remaining after n iterations of the Farey map. A similar exponential decay can be seen for 
the full dynamics if we use the T time coordinate. Naturally, the rate of decay is coordinate 
dependent in general relativity. However, <r> is related to the Lyapunov dimension which 
is coordinate independent. 

To build the pockets we modify the map to, 



oo, 



(exit) 



__ J u n - 1, 2 < u n < u cxit (oscillation) 

1 



(3.15) 



u„. — l' 



1 < u n < 2 



(bounce) . 



We refer to this map as the Farey exit map. 

In Fig. 8, we show a histogram of the fixed points of the exit map ( p.!5|) by finding all 
roots up to order k — 16 with u cx it = 8. As promised, the strange repellor for the Farey exit 
map is little changed from the repellor of the full map shown in Fig. 1. The main difference 
is that the gaps around the rationals are now completely devoid of roots and are not just 
sparsely populated regions. As a result, the capacity dimension of the fractal in Fig. 8 will 
not saturate at the phase space dimension. However, since the dense regions of the repellor 
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have not been affected by the introduction of exits we expect the information dimension to 
be little changed from that of the full map. We find = 1.87 ± 0.01, in agreement with 
the information dimension of the future invariant set shown in Fig. 2. 



C. Measures and <A> 

When dealing with maps, the transient lifetime and Lyapunov exponents for the aperiodic 
trajectories can be calculated from a knowledge of the periodic trajectories. In this way, the 
strange repellor provides a complete description of the dynamics even though the periodic 
trajectories are a set of measure zero. 

Recall that the fixed points Ujk at order k are found by solving the 2 k — 1 equations 
F k (ujk) = Ujk- Writing the stability coefficient of Uj k as Cjk, it can be shown that 

H(S)= lim £ — , (3.16) 

where the measure covers the region S. If S covers the entire invariant set we find 

1 ( k \ 
fi = lim V" — = lim exp . (3-17) 

fc-*oo j Cjk k^oo \ <T> J 

If <t> 7^ oo, the measure decays as trajectories are lost from the system. The non-negative 
Lyapunov exponent is given by |39| 



<A>= lim E > C ? )/c '- fc . (3.18) 

To clarify, < A > is the average Lyapunov exponent for a typical aperiodic trajectory. This 
is to be contrasted with the average Lyapunov exponent, A p , of the rare, periodic orbits 
described in § [11 C\ 

Applying these relations to the Farey map without exits we find 

fi k = — - 1 - e- a5 fc- a268 , (3.19) 

so that <t> = oo. As emphasized earlier, the Farey map without exits sees typical orbits 
return to the scattering region so \i = 1. 

Before we calculate D^, we take a brief aside to demonstrate the effectiveness of the 
measure (|3.16|) by reproducing the topological entropy calculated analytically in Eqn. ( |2.21|) . 



Since the measure does not decay, we can use the g-order entropy spectrum [K| . 

y 3 

as an independent calculation of the topological entropy. The topological entropy is recov- 
ered in the limit q — > 0, so that H T = H . Using all roots up to k = 16 and taking q = 0.01 
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we find i/o.oi — 0.693 ~ In 2, in good agreement our analytic calculation of Ht- The fast 
convergence of Hqqi as a function of k is shown in Fig. 9. In theory, by taking the limit 
q — > 1 we should recover the metric entropy h(fi) = Hi. In practice, the expression for Hi 
converges very slowly with k and we are unable to confirm that h = Hi = 0. 

Turning to the Farey exit map with an exit at w ex it = 8 we find from Eqns. Q3.17D and 
( gig ) that <r> = 12.6 ± 0.1 and <A> = 0.724 ± 0.005. The convergence as k -* oo of 
Eqns. ( |3.17| ) and ( |3. are shown graphically in Figs. 10 and 11 respectively. Notice that 
the fits only start to converge once k > u exit + 1. This makes sense as the the largest root at 
order k is w max — k — 1. The map is unaware it has exits until orbits start to escape. These 
results for <r> and <A> combine to tell us that D L = 1.78 ± 0.02 and h = 0.645 ± 0.005. 

As an independent check we can use a direct Monte-Carlo evaluation based on evolving 
a random collection of initial conditions. The Monte-Carlo method allows us to estimate 
the properties of typical trajectories by measuring the properties of a large random sample 
then infering from this the properties of the collective. 

By evolving iV(0) = 10 5 randomly chosen initial points, we find from Eqns. (|3.2|) and 
(0) that <r> = 12.2 ±0.1 and < A> = 0.74 ±0.01. The convergence of <r> and < A> as 
a function of the number of iterations n is shown in Figs. 12 and 13. In each case the initial 
conditions in the first run were taken from the interval u = [1, 2], and those in the second 
run from the larger interval u — [1,8]. In this way, the asymptotic value was approached 
from above and below, leading to a more accurate estimate of the limit. After about n = 60 
iterations statistical errors start to become large since the number of points remaining in 
the interval becomes small. This is apparent in Figs. 12 and 13. Using the Monte-Carlo 
values for < r > and < A > we find Dl = 1.78 ± 0.03 and h = 0.66 ± 0.01. The two methods 
agree. 

In order to compare Dl to the information dimension of the strange repellor we need to 
numerically generate the repellor and find its fractal dimensions. Recall that in Fig. 8 we 
generated the unstable manifold of the Farey exit map by solving for all roots up to order 
k = 16. With an exit at M ex it = 8 we find D\ = 1.87 ± 0.01. Consequently, the information 
dimension of the strange repellor is Di = 1.74 ± 0.02, which agrees with the result without 
exits (Eqn. ( |2.32j )). As mentioned previously, the dimension found for k = 16 is a lower 
bound for the true dimension. While there is reasonable agreement between Di and Dl 
at order k = 16, we would like to test the relation at higher orders in k. However, this is 
difficult since the number of roots increases exponentially with k. 

Alternatively, we can expose the future invariant set by the Monte-Carlo method. By 
starting off N(0) = 5 x 10 6 points in the interval u = [1, 2] and iterating them n = 60 
times, we are left with N ~ 3 x 10 4 points. These remaining points must closely shadow 
the future invariant set. A histogram of the set generated in this way is shown in Fig. 14. 
The gap structure is identical to that seen in Fig. 8 using all roots at order k = 16, but the 
bins are more evenly filled. As a result, the information dimension takes the slightly higher 
value -D" = 1.89 ±0.01. This method gives an information dimension of the strange repellor 
closer to Di = 1.78 ± 0.02. Since the Monte-Carlo method only finds orbits that shadow 
the repellor, and not the true periodic orbits, we expect the fractal dimension found in this 
way will be an upper bound. 

As expected, the results for < r > and < A > depend on how wide we make the pocket. 
In the limit u exit — > oo we know <r >— > oo and < A>— > 0. It would be interesting to see if 
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the product <r><A> remains finite as w ex it becomes large. Unfortunately, the numerical 
values converge very slowly as the pockets become thinner so we were unable to test whether 
D L remains fixed as M exit — > oo. 

To conclude this section, we find from two separate methods that 

D L = 1.78 ±0.03 . (3.21) 

Calculating the information dimension in two different ways gives the range 

Di = (1.74 -> 1.78) ± 0.02 . (3.22) 

Systematic errors and slow convergence of the dimension can account for the discrepancy 
in the two values for D\. The Lyapunov dimension is in reasonable agreement with two 
different calculations of the information dimension. The Lyapunov dimension thus fares 
well as a coordinate invariant combination of the otherwise ambiguous Lyapunov exponent 
and metric entropy. 



IV. COMPARING TO THE GAUSS MAP 

Barrow used the Gauss map developed by BKL to derive dynamical systems properties 
such as metric and topological entropies. In this section we compare Barrow's results for 
the Gauss map to our results for the Farey map. 

There are two critical features which distinguish the treatment of the maps. Firstly, the 
Farey and Gauss maps are topologically inequivalent. Secondly, the techniques used to study 
the maps are different. We begin by isolating the strange repellor before using the repellor 
to reconstruct the properties of typical aperiodic trajectories. In contrast, Barrow worked 
directly with the aperiodic trajectories of the Gauss map as opposed to the complement, 
the repellor. 

The Gauss map can be obtained from the Farey map by chopping out the oscillations. 
This was first done by BKL when they realised that the sensitive dependence on initial 
conditions was due to the bounces and not the oscillations. The Gauss map emerges from 
the Farey w-map as follows: If u n > 2 then oscillations will take place until u n < 2. The 
bounce sequence that follows can then be expressed as 

1 . . 

un+i = f — 7 (4.1) 

U N - [Un\ 

where N is the number of bounces. Eqn ( |4.1| ) is just the Gauss map, G(u). The Gauss map 
is not invertible since it is derived from F(u) and not the full F(u, v). 

The Farey and Gauss maps correspond to different discretisations of the full dynamics. 
The discrete ticks of the Farey clock correspond to periodic behaviour in terms of the con- 
tinuous time variable T = lnln(l/t). This is not true for the Gauss map. The time between 
bounces for the Gauss map can vary dramatically relative to T time. In compactifying the 
discrete time coordinate, tp = n, of the Farey map to arrive at the discrete time coordi- 
nate, tc = N, of the Gauss map, holes have had to be cut. There is no smooth coordinate 
transformation connecting the time variables tp and to- As a result, the F and G maps are 
topologically inequivalent and have different topological entropies. 
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A. The repellor 



Although topologically inequivalent, there are deep similarities between the two maps F 
and G. Fore one, they have identical fixed points. Still, some features of the repellor are 
altered. 

Consider orbits with a Gauss period Pq = 1, so that the fixed points are obtained from 



The solutions are 



1 



Ui 



Ux - Mi 



(4.2) 



Ux 



nix + Jmf + 4 



(4.3) 



where again 



mi 



(4.4) 



Notice, these ux correspond to the silver means of Eqn. ( |2.15 ). Hence, the period Pq = 1 
orbit of the Gauss map generates all of the single bounce orbits of the Farey map, but with 
the longer period Pp = rax- Notice also that there are already an infinite number of fixed 
points for the Gauss map at order N — 1. 

We can calculate the Lyapunov exponents for the Pq = 1 orbits: 



\ G (u)= lim Tr Hu) 2N 



2 In 



mx + yrrix + 4 



(4.5) 
(4.6) 



And so the silver mean exponents for the two maps are related by 

\ G (u) = mx \ F (u) . 



(4.7) 



The exponent for G is larger than that for F. 

Continuing in the same vein, consider the period Pq = 2 fixed points. The condition 



G 2 (u) 



generates the doubly infinite set of roots 



u . 



(4.8) 



u 



mx + J m\ + 4mi/m 2 
2 

2 



(4.9) 



m\ + 4mi / 777,2 — TTl\ 



Notice that all the P& = 1 orbits are included in this solution when mx = m 2 . By considering 
the Gauss map we have been able to solve for all of the two-bounce fixed points of the Farey 
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map at order k = m\ + ni2. These periodic orbits have 2 bounce epochs separated by an 
m\ — 1 oscillation era and an m 2 — 1 oscillation era. The Lyapunov exponents for the period 
2 orbits are given by 



\ G (u)=\n\ m ° +1[ . (1.10) 

'm\ + Ami/rri2 

These are related to the Lyapunov exponents for a corresponding periodic orbit of the Farey 
map by 

AO(fl) = ^ + n^ m (4 n) 

In terms of the cfe, the period of the Gauss map is the length of the generator for the 
periodic irrational. For example, the periodic irrational [1, {mi, 1712, . . . , m^}] has Pq = N 
and Pf = J2i=i m %- From this we can deduce the general result 

X G _P F _T^m i 

X F ~P G ~ N • [4AZ) 

Interestingly enough, while the Lyapunov exponents and the period of the two maps differ, 
their product is invariant: 

X G P G = X F P F . (4.13) 

Since typical periodic orbits of the Farey map have equal numbers of bounces and oscil- 
lations, it follows that the average Lyapunov exponents for the periodic orbits satisfy 

A G = 2\ F ~ 1.59. (4.14) 

A direct calculation of A G confirms this expectation. 



B. Comparing entropies 

In terms of continued fractions it is a simple matter to write down all the periodic orbits 
of the Gauss map. The period 1 orbits are given by u = [{mi}], the period 2 orbits by 
u = [{mi, m2}] etc. From this it follow that the number of orbits with period Pq = K scales 
as N(K) = oo K . Because the number of roots at each order is infinite, a regularisation 
procedure has to be introduced before quantities such as H T and h(fM) can be calculated. 

After employing a suitable regulator M, the measure ( |3.16| ) can be expressed in terms 
of an integral over x = 1/u: 



s 



li(S) = I p(x)dx, (4.15) 

where 
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Taking S to cover the entire region x = [0, 1] (u = [1, oo]), we have fj, = 1 and < r > = oo. 
From the Gauss map in terms of the x coordinate, G(x) = x^ 1 — [x" 1 ], we find G'(x) = —x~ 2 
and 



rl rl — 2 hi £ 

<A> = / In \G'(x)\p(x)dx — / — ; — dx 

Jo 1 V JlHK J Jo (1 + x ln2 



7T 2 

61n2 



(4.17) 



This expression is smaller by a factor of In 2 from the usual expression. The discrepancy 
can be traced to the use of log 2 rather than natural logarithms in Ref. H. Similarly, the 
topological entropy of the Gauss map is given by 

HS = -f-, (4.18) 
T 61n2' V ; 

where again our expression differs from the oft quoted result by a factor of In 2. As promised, 
the topological entropy of the Gauss map greatly exceeds that of the Farey map. Moreover, 
the metric entropy of the Gauss map, h G =<\> — <t>~ 1 = 7r 2 /(61n2) is finite while the 
metric entropy of the Farey map vanishes. 



C. Fractal dimension of the Gauss repellor 

Although the Gauss map and the Farey map share common fixed points, they differ in 
how the fixed points are distributed as a function of period. It is impossible to produce a 
histogram in analogy to that of Fig. 1 since the Gauss map generates an infinite number of 
fixed points at every order. We can however argue that its very likely the fractal generated 
by the Gauss map is uniform and not multifractal. In other words, all the D q 's are likely the 
same. Indeed, we can argue that the fractal dimensions saturate the phase space dimension, 
D q — D — 1, for all q. The reason being that the periodic irrationals are dense on the 
number line. We have already argued in §[11 A| that this ensures D a = 2 for the Farey map. 



Let us consider this argument more closely. As Pq becomes large, the number of terms 
in the cfe cycle becomes large also. The longer the cfe cycle, the more uniform the coverage 
of the number line. For orbits of period P, the gaps between roots are always less than 
1/ fp, where 

f P = ^= ( 7 ^) + (-1) V (P+1) ) • (4-19) 



v/5 



Here 7 = (v^5 + l)/2 is the golden mean. For example, by the time P = 10, fio = 89 and 
the largest possible gap is bounded by ~ 0.0001. In comparison, the Farey map at order 
k = 10 has a finite number of roots and gaps as large as ~ 0.2. Moreover, the gaps close 
exponential fast with P for the Gauss map, while they close slower than 1/k for the Farey 
map. While the rate of closure makes no difference to D a , the multifractal dimensions with 
higher g's take into account the density of points. The Gauss map lays down a very even 
covering of the numberline, so we expect D Q = D\ etc. 
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If the conjecture that D\ = Dl holds true, then we should find Dl = 1. This is in fact 
what we find. For the 1-D map 



D L = \ . (4.20) 

<A><r> 

The Gauss map is characterized by <r> = oo and <A> = 7r 2 /(61n2). Inserting these values 
into ( 4.20|) we find Dl = 1, in accordance with our argument that Di — 1. 



As with the topological entropies, we find the fractal dimensions of the two maps are 
inequivalent. This does not invalidate the statement that fractal dimensions are coordinate 
invariant signals of chaos. The two maps cannot be connected by a smooth coordinate trans- 
formation. They are, as we have already argued, topologically inequivalent. The discrete 
Farey time is topologically the same as the T = lnln(l/t) of the full dynamics. The discrete 
Gauss time on the other hand changes the topology of time by chopping out the oscillations. 



V. NUMERICAL RESULTS FOR THE FULL DYNAMICS 

The Farey map has taught us a lot about the mixmaster universe. To be certain the 
Farey map provides an accurate description, we now expose the strange repellor in the 
minisuperspace of the full continuum dynamics. 

Continuous dynamical systems are much harder to solve than discrete dynamical systems. 
We cannot hope to find the periodic orbits analytically. Instead, we apply the same kind of 
Monte Carlo approach used to check the analytic results for the Farey map. Our ambition 
will be limited here to finding the fractal manifest in minisuperspace phase space and the 
fractal dimension. 

Exposing the chaotic invariant set in a dissipative system is simple. By randomly choos- 
ing a collection of initial conditions and evolving them according to the dynamical equations, 
the trajectories will soon settle onto the various attractors in phase space. The decrease in 
available phase space volume forces all trajectories onto the chaotic attractors. For a Hamil- 
tonian system this is not the case. Phase space volume is conserved and the chaotic invariant 
set remains hidden among the slurry of aperiodic trajectories. Nonetheless, the chaotic in- 
variant set can be uncovered by a variety of methods. For compact Hamiltonian systems the 
proper interior maxima (PIM) procedure MM efficiently reveals the strange repellor while for 



non-compact Hamiltonian systems, including Hamiltonian exit systems ||41|1 , fractal basin 
boundaries are the prefered method. 



A. The PIM procedure 

The proper interior maxima (PIM) procedure is able to isolate the strange repellor in 
most Hamiltonian systems. In particular, it could be applied to the mixmaster system with 
or without exits. Although we did not apply the PIM procedure to the mixmaster, we 
mention it here as an alternative to introducing exits. 

The basic idea behind the PIM procedure is ingeniously simple. As discussed earlier, the 
stable and unstable manifolds of the strange repellor are interchanged under time reversal. 



29 



As time moves forward trajectories near the repellor are repelled along the unstable man- 
ifold and attracted along the stable manifold. Evolving the system into the past reverses 
the attracting and repelling directions. The PIM procedure exploits this property in the 
following way: a) Evolve a collection of trajectories near the strange repellor into the future 
and into the past by an equal amount of time, b) Isolate subsets of the evolved bundles that 
remain nearby in both the past and the future and discard the remainder, c) Zero in on 
the surviving trajectories and form a new, smaller bundle of trajectories around them and 
repeat the procedure. After a few iterations the surviving trajectories will closely shadow 
the invariant set as noninvariant trajectories have been discarded. 

The PIM procedure is shown schematically in Fig. 15. To simplify the picture, the initial 
bundle of trajectories is shown centered on a point belonging to the strange repellor (black 
dot). The labels a) and b) correspond to the steps in the PIM procedure. The great feature 
of the PIM procedure is the way it uses the instability of the repellor's orbits to its advantage. 
Consider the original set of trajectories in the circle of radius R about the invariant point. 
Evolving into the future and into the past by an amount At leads to ellipses with axes of 
length i?exp(— <A> At) and i?exp(+ <A> At). These ellipses share a At invariant set 
of radius r ~ i?exp(— <A> At). Thus, after a few iterations, the repellor lies exposed. In 
effect, the PIM procedure turns the strange repellor into a strange attractor. 

Using a variant of the PIM procedure entire trajectories belonging to the strange repellor 
can be reconstructed. While it would be valuable to apply these techniques to the mixmaster 
system, the numerical implementation of the procedure is difficult. A much simpler method 
of exposing the strange repellor is to introduce exits, and then chart the fractal coastline of 
the outcome basins. 

B. Fractal basin boundaries 

Rather than hunting for the chaotic skeleton among the aperiodic trajectories, the Hamil- 
tonian exit method creates a sieve through which the slurry of aperiodic trajectories escape, 
leaving the strange repellor in clear view. Like the PIM procedure, the introduction of exits 
exploits the instability of the chaotic trajectories. The faster trajectories are expelled by the 
repellor, the faster the invariant set is exposed. 

The repellor is manifest at the boundary in phase space which separates initial conditions 
on the basis of their outcomes. These basin boundaries often become fractal in chaotic 
scattering. The fractal basin boundary is composed of trajectories which spend substantial 
time on the repellor. 

For the full mixmaster dynamics we employ the lessons of the Farey map. Indeed, it is a 
simple matter to directly implement exactly the same exit conditions we used for the map. 
By considering the ratios (lna)'/(ln&)' et eye. (a,b,c), we can test to see if a mixmaster 
universe is coasting in a Kasner phase. If it is, we can then use Eqn. ( |1.5| ) to read off the 
value of u. When u > u ex a we terminate the evolution and assign an outcome based on which 
axis is collapsing most quickly. By colour coding the initial conditions near the maximum 
of expansion according to their outcome near the big crunch singularity, we produce plots 
of the outcome basins. 

Using a fourth order Runge-Kutta integrator with adaptive step size, we evolved 300 x 
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300 grids of initial conditions and recorded their outcomes. The Hamiltonian constraint, 
Eqn. (|1.1|) was monitored at all times, and the error tolerances were adaptively corrected to 
ensure the constraint was satisfied to within 1 part in 10 5 along each trajectory^. In order 
to make contact with the Farey map, the initial conditions were chosen by fixing uj = 1/3 
and taking (u a , v Q ) from a 300 x 300 grid. The value of Q Q was then found by solving the 
Hamiltonian constraint. Inserting these values of (u Q , v Q , Q a , uj ) into Eqn. (|1.3j ) yielded the 
initial conditions for the numerical integration of the equations of motion. 

A portion of the basin boundaries in the (u, v) plane is displayed in Fig. 16. Depending 
on which axis is collapsing most quickly when the trajectory escapes, the initial grid point 
is coloured black for a, grey for b and white for c. The exit was set at w exit = 8 in order 
to make comparison with the Farey exit map results. The basins form an intricately woven 
tapestry of roughly vertical threads. The basin boundaries appear to form a cantor set of 
vertical lines. By zooming in on a small portion of Fig. 16, we see from Fig. 17 that the 
dense weave persists on finer and finer scales. 

Points belonging to the fractal basin boundaries comprise a future invariant set. Tra- 
jectories belonging to this set never decide on a particular outcome, and so never escape 
through an exit. These are the trajectories that form the unstable manifold of the strange 
repellor. It is instructive to compare the future invariant set seen in Figs. 16 and 17 with 
the Farey map's future invariant set shown in Fig. 2. The pattern of gaps and dense regions 
is strikingly similar. The warpage of the vertical stripes in Figs. 16 and 17 can be accounted 
for by our choice of starting point. Because our initial conditions start the universe near the 
maximum of its expansion, we are looking at the future invariant set in a region where the 
approximations used to derive the Farey map are very poor. Nonetheless, as the trajectories 
evolve toward the big crunch the continuum dynamics settles onto a pattern well described 
by the Farey map. Since the fine structure of the weave is laid down after many bounces 
and oscillations, the fractal structure produced by the full equations should be much the 
same as that produced by the Farey map. 

In order to test this proposition we numerically evaluate the fractal dimension of Fig. 17. 
To do this we employ the uncertainty exponent method described in §11 A, Randomly 



choosing 1000 points in the region covered by Fig. 17, we record how many have certain 
outcomes as a function of initial uncertainty 5. A plot of \nf(5) as a function of In 5 is 
shown in Fig. 18. According to Eqn. ( |2.30| ), the gradient of this graph yields the uncertainty 
exponent a = 0.14 ± 0.01 and a capacity dimension of -D" = 1.86 ± 0.01. 

The method described above delivers the capacity dimension of a sample of the repellor. 
The capacity dimension of a sample often is not equal to the capacity dimension of the full 
repellor. Rather, the capacity dimension of the sample equals the information dimension of 
the full repellor |4"5 |. The reason is the following: The random sampling of the Monte-Carlo 
approach favours the densest regions of the fractal. The method naturally weights therefore 
the densest regions. Since the weighted dimension is in fact the information dimension, 
it follows that of the sample actually equals D\. We have therefore really found the 
information dimension of the mixmaster's future invariant set: 



4 We check the normalised constraint H = H/ [o 4 + 6 4 + c 4 + 2(b 2 c 2 + a 2 c 2 + b 2 a 2 )] . 
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D\ = 1.86 ± 0.01 



(5.1) 



The information dimension of the repellor from these numerical experiments is then given 
by 

D 1 = 2D'1 - D = 1.72 ±0.02 . (5.2) 

Within numerical uncertainty, this number agrees with what we found for the Farey map 
(Eqns. ( p.32j ) and ( |3.22| )). The small discrepancy is probably due to the approximations 



used in deriving the map, or systematic numerical errors. 

It is worth noting that the fractal basin boundaries can be uncovered in any 2-D slice 
through the 6-D phase space. For example, in Fig. 19 we display the outcome basins in 
the anisotropy plane (/3+,/3_). The initial conditions were chosen by setting Q = —1, a = 
10, $ = 10, and selecting (/3 + ,/3_) from a 600 x 600 grid. The remaining coordinate, 7, 
was fixed by the Hamiltonian constraint. Again we see an interesting mixture of regular 
and fractal boundaries, this time forming a symmetrical patchwork about a central axis of 
symmetry. The symmetry of the basins reflects the symmetry of the MSS potential shown 
in Figs. 6 & 7. 

By choosing different slices or different coordinate systems, many different views of the 
future invariant set can be uncovered. However, these are purely cosmetic changes. No 
matter what coordinates we choose there will always be fractal basin boundaries, and these 
boundaries will always have the same fractal dimension. No matter how you look at it, the 
mixmaster universe is chaotic. 



VI. SUMMARY OF RESULTS, DISCUSSION 

The power of the fractal in relativity is its observer independence. We isolate two related 
fractals. The strange repellor of the Farey map is shown to be a multifractal and analyzed 
in detail. The fractal repellor is then excavated numerically in the phase space of the full, 
unapproximated dynamics. 
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The full force of the conclusions comes with a comparison of the dimensions found in 
entirely different manners. We have collected our results for the discrete time maps in the 
chart above. For comparison with the chart, the information dimension of the fractal basin 
boundary in the full dynamics is D\ = 1.72 ± 0.02. Not only is a fractal a coordinate 
invariant declaration of chaos, we have further found that the three information dimensions 
agree within errors: the information dimension of the repellor from the discrete Farey map, 
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the Lyapunov dimension of the Farey map, and the information dimension of the fractal 
basin boundaries. As has long been suspected, the collective of mixmaster universes is 
certainly chaotic. 

In a sense, chaos in the mixmaster universe makes it impossible to create one in the first 
place. The mixmaster is anisotropic but homogeneous. The extreme sensitivity to initial 
conditions induces an instability to inhomogeneities. Two different regions of spacetime 
with even slightly different initial conditions quickly evolve away from each other [J46|] . The 
universe then becomes inhomogeneous as well as anisotropic, rendering the big crunch similar 
to a generic inhomogeneous collapse to a black hole. The churning dynamics on approach 
to a singularity may therefore be temporally and spatially chaotic, whether it be in a dying 
star or a dying cosmos. 

Our gratitude to J. Barrow, C. Dettmann, N. Frankel, P. Ferreira and J. Silk for their 
comments and valuable discussions. We also thank C. Dettmann for letting us adapt his 
computer codes. N.C thanks Case Western Reserve University for their hospitality and 
generous use of their computers. 
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FIGURES 
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FIG. 3. Finding the information dimension of the future invariant set. 
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FIG. 4. The strange repeller in the region 1 < u < 2, 1 < v < 2. 
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FIG. 6. Equipotentials of the minisuperspace potential. 
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FIG. 9. Numerical convergence of the topological entropy, -ffo.oi — Ht, as a function of k using 
Eqn. (3.20). The solid line indicates log 2. 
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FIG. 10. Finding the decay time <r> for the Farey exit map using Eqn. (3.17). The straight 
line fit yields a decay constant of <r>= 12.6 ± 0.1 
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FIG. 11. Finding the Lyapunov exponent <A> for the Farey exit map using Eqn. (3.18). The 
straight line fit yields a Lyapunov exponent of < A>= 0.724 ± 0.005. 
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FIG. 12. Finding <A> for the Farey exit map using the Monte Carlo method. 
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FIG. 13. Finding <r> for the Farey exit map using the Monte Carlo method. The decay time 
is found from the slope of the line of best fit. 
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FIG. 14. A histogram of points shadowing the Farey exit map's future invariant set found by 
Monte Carlo methods (u cx it = 8). 
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FIG. 17. A portion of Fig. 16 magnified 50 times. 
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FIG. 18. Calculating the uncertainty exponent for Fig. 17. 
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FIG. 19. The outcome basins in the anisotropy plane. 



45 



